## R code for a behavioral model of turnout
## James Fowler
## December 8, 2003
##
## Note: This is the stripped down code.  If you want any data from
## this you'll need to add some lines to print to the screen or a file.
##

rm(list=ls(all=TRUE))

nPeriods<-1000 ## number of periods
nSims<-1000 ## number of simulations

nDems<-5000 ## number of Democrats
nReps<-5000 ## number of Republicans
n<-nDems+nReps ## number of citizens

winPayoffD<-1.0 ## Dem payoff for winning
winPayoffR<-1.0 ## Rep payoff for winning

losePayoffD<-0 ## Dem payoff for winning
losePayoffR<-0 ## Rep payoff for winning

costD<- 0.25  ## cost to the democrats
costR<- 0.25 ## cost to the republicans

iaspirationD<- 0.5 ## initial aspiration Dems
iaspirationR<- 0.5 ## initial aspiration Reps

iturnoutpropensityD<-0.5 ## initial propensity to turnout Dems
iturnoutpropensityR<-0.5 ## initial propensity to turnout Reps

## auxilary parameters

tau<- 0 ## if 1, use Bush Mosteller rule, if 0, unbiased adjustment

alpha<-0.1  ## propensity update weight for success
beta<-0.1   ## propensity update weight for failure
lambda<-0.95 ## weight for aspiration update
inert<-0.01 ## probability a voter updates propensity or aspiration
support<-0.2 ## support of random payoff shock

## payoff function
payofff<-function(winner,preference,cost,action) 
  preference*(winner*winPayoffR+(1-winner)*losePayoffR)+
  (1-preference)*(winner*losePayoffD+(1-winner)*winPayoffD)-
  action*cost+round(runif(length(action),-support/2,support/2),3)

## aspiration update function
aspirationf<-function(aspiration,payoff)
  ((aspiration>payoff)*floor(1000*(lambda*aspiration+(1-lambda)*payoff))+
  (aspiration<payoff)*ceiling(1000*(lambda*aspiration+(1-lambda)*payoff)))/1000+
  (aspiration==payoff)*aspiration

## adjustment for propensities
  propensityf<-function(propensity,aspiration,action,payoff)
    pmin(1,pmax(0,((action) *
     ((payoff>=aspiration)*ceiling(1000*(propensity+alpha*(1-tau*propensity)))+
      (payoff<aspiration)*floor(1000*(propensity-beta*(1-tau*(1-propensity)))))+
    (1-action) * 
     ((payoff>=aspiration)*floor(1000*(propensity-alpha*(1-tau*(1-propensity))))+
      (payoff<aspiration)*ceiling(1000*(propensity+beta*(1-tau*propensity)))))/1000))

## voter preferences
preferences<-c(rep(0,nDems),rep(1,nReps))

## voter costs
costs<-c(rep(costD,nDems),rep(costR,nReps))

for (i in 1:nSims) {

  ## initial values
  voteDts<-NULL ## init Dem vote time series
  voteRts<-NULL ## init Dem vote time series

  ## propensity to vote
  propensities<-c(rep(iturnoutpropensityD,nDems),rep(iturnoutpropensityR,nReps))

  ## voter aspirations
  aspirations<-c(rep(iaspirationD,nDems),rep(iaspirationR,nReps))

  for (j in 1:nPeriods) {

    ## voters choose whether or not to vote
    actions<-runif(n,0,1)<propensities

    ## total votes
    voteD<-sum((1-preferences)*actions)
    voteR<-sum(preferences*actions)
    voteDts<-c(voteDts,voteD/nDems)
    voteRts<-c(voteRts,voteR/nReps)

    ## winner determined
    winner<-(voteR>voteD)+(voteD==voteR)*(runif(1,0,1)>0.5)

    ## voter payoffs
    payoffs<-payofff(winner,preferences,costs,actions)

    ## set up voter index vector for propensities
    noninertials<-which(runif(n,0,1)>inert)

    ## update propensities
    propensities[noninertials]<-propensityf(propensities[noninertials],
      aspirations[noninertials],actions[noninertials],payoffs[noninertials])

    ## set up voter index vector for aspirations
    noninertials<-which(runif(n,0,1)>inert)

    ## update aspirations
    aspirations[noninertials]<-aspirationf(aspirations[noninertials],
      payoffs[noninertials])
  }
}